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ABSTRACT 

Considerable progress has been made in recent years in observations of atmospheric signa¬ 
tures of giant exoplanets, but processes in rocky exoplanets remain largely unknown due to 
major challenges in observing small planets. Numerous efforts to observe spectra of super- 
Earths, exoplanets with masses of 1-10 Earth masses, have thus far revealed only featureless 
spectra. In this paper we report a 4-a detection of variability in the dayside thermal emis¬ 
sion from the transiting super-Earth 55 Cancri e. Dedicated space-based monitoring of the 
planet in the mid-infrared over eight eclipses revealed the thermal emission from its day- 
side atmosphere varying by a factor 3.7 between 2012 and 2013. The amplitude and trend 
of the variability are not explained by potential influence of star spots or by local thermal or 
compositional changes in the atmosphere over the short span of the observations. The possi¬ 
bility of large scale surface activity due to strong tidal interactions possibly similar to Io, or 
the presence of circumstellar/circumplanetary material appear plausible and motivate future 
long-term monitoring of the planet. 
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1 INTRODUCTION 

A key question in exoplanet science is about the origin of planets in 
tight orbit (P<0.75 days) around their stars. During its primary mis¬ 
sion, Kepler enabled the detection of a hundred ultra-short period 
(USP) planet candidates (Sanchis-Ojeda et al. 2014). The majority 
of these USP planets are small (Rp < 27?®) and may suggest that 
they have experienced a phase of intense erosion due to tidal activ¬ 
ity and/or intense stellar irradiation (e.g. Owen & Wu 2013). 

One example of these disintegrating USP planets is 
KIC12557548b (Rappaport et al. 2012), which exhibits changes 
in the optical transit depth and shape that seem consistent 
with an evolving cometary-like environment (Brogi et al. 2012; 
Rappaport et al. 2014; Croll et al. 2014; Bochinski et al. 2015). 
Scattering of dust particles surrounding the planet are thought 
to be the driver of the observed variability, as the engulfed sub- 
Mercury-size object itself would be too small to be actually de¬ 
tected (Perez-Becker & Chiang 2013). The high effective tempera¬ 
tures of USP planets could enable a silicate-rich atmosphere made 
of refractory material to be sustained, which would shape the over¬ 
all opacity. Replenishment of grains in the atmosphere could be 
achieved by condensation at cooler temperatures near the day- 
night terminator or at lower pressures (Schaefer & Fegley 2009; 
Castan & Menou 2011). Another source for the observed refractory 
material could be volcanism, fuelled by tidal heating (Jackson et al. 
2008; Barnes et al. 2010) similar to Io, where the interactions be- 
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tween an inhomogeneous circumplanetary neutral cloud and a cir- 
cumstellar charged torus (Belcher 1987) would result in transit 
depth variability over timescales that are similar to the planet’s or¬ 
bital period. 

The super-Earth 55 Cancri e has a mass of 8.09 ± 0.26 M® 
(Nelson et al. 2014), a published radius of 2.17 ± 0.10 7?® 
(Gillon et al. 2012), and orbits a nearby sun-like star with an 18- 
hour period. 55 Cnc e is one of the largest members of known USP 
planets. Owing to its close proximity to the parent star, 55 Cnc e 
has an equilibrium temperature of ~2400 K and may experience 
intense tidal heating (Bolmont et al. 2013). With V~6, 55 Cnc 
is the brightest star known to host a transiting exoplanet, and 
has led to measurements of the planet’s radius at exquisite preci¬ 
sion in the visible as well as in the Spitzer 4.5 pm IRAC band 
(Demory et al. 2011; Winn et al. 2011; Gillon et al. 2012). Ther¬ 
mal emission measurements in the Spitzer 4.5 /im IRAC photo¬ 
metric band (Demory et al. 2012) yielded a brightness temperature 
of about 2,300 K in this channel. Such elevated temperature has 
implications on the nature of the planet atmosphere. First, this high 
temperature means that for a given chemical composition, 55 Cnc e 
would have a larger atmospheric scale height compared to a cooler 
planet with the same composition. Second, clouds are less likely to 
form in the atmospheres of extremely hot exoplanets compared to 
cooler objects (Madhusudhan & Redfield 2015), which is particu¬ 
larly relevant in the light of the numerous non-detections of spec¬ 
tral features in exoplanets (e.g. Pont et al. 2008; Sing et al. 2011; 
Knutson et al. 2014; Kreidberg et al. 2014). The combination of the 
high photometric precision possible for transit/occultation obser- 
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vations and the high planetary temperature means that 55 Cnc e 
currently represents one of the most promising super-Earth for in¬ 
frared spectroscopic observations of its atmosphere. Using the pub¬ 
lished mass and radius alone, the interior composition of the planet 
is consistent with: (1) a silicate-rich interior with a dense H 2 O en¬ 
velope of 20% by mass (Demory et al. 201 1), (2) a purely silicate 
planet with no envelope (e.g. Gillon et al. 2012), or (3) a carbon- 
rich planet with no envelope (Madhusudhan et al. 2012) among 
other compositions. 

In this paper, we present an analysis of archival and new 
Spitzer/ IRAC data obtained in 2012 and 2013. This combined 
dataset suggests a radically different nature for 55 Cnc e that could 
be more consistent with a scenario similar to KIC 12557548b than 
the more conventional picture that was attached to 55 Cnc e and 
which is described above. We present the observations and data 
analysis in Section 2, a comparison with previous results in Sec¬ 
tion 3 and we finally propose an interpretation of our results in 
Section 4. 

2 DATA ANALYSIS 

2.1 Observations 

This study is based on six transits and eight occultations of the 
super-Earth 55Cnee acquired between 2011 and 2013 in the 4.5- 
//m channel of the Spitzer Space Telescope (Werner et al. 2004) 
Infrared Array Camera (IRAC) (Fazio et al. 2004). Data were ob¬ 
tained as part of the programs IDs 60027, 80231 and 90208. 
Two transits were observed in 2011 (PID 60027) (Demory et al. 
2011; Gillon et al. 2012), four occultations in 2012 (PID 80231) 
(Demory et al. 2012), as well as four transits and four other oc¬ 
cultations in 2013 (PID 90208). The corresponding data can be 
accessed using the Spitzer Heritage Archive 1 . All Spitzer obser¬ 
vations have been obtained with the same exposure time (0.02 s) 
and observing strategy (stare mode). The last three occultations of 
PID 80231 have been obtained using the PCRS peak-up mode as 
all of our 2013 data. The PCRS peak-up procedure is designed to 
place the target on a precise location on the detector to mitigate the 
intra-pixel sensitivity (Ingalls et al. 2012, and references therein). 
The AORs 48072704 and 48072960 experienced 30-min long inter¬ 
ruptions that happened during the ingress and egress of 55Cnee’s 
transit respectively. 

2.2 Data reduction 

The starting point of our reduction consists of the basic calibrated 
data (BCD) that are FITS data cubes encompassing 64 frames of 
32x32 pixels each. Our code sequentially reads each frame, con¬ 
verts fluxes from the Spitzer units of specific intensity (MJy/sr) to 
photon counts, and transforms the data timestamps from BJDutc 
to BJDtdb (Eastman et al. 2010). The position on the detector of 
the point response function (PRF) is determined by fitting a Gaus¬ 
sian to the marginal X, Y distributions using the GCNTRD pro¬ 
cedure part of the IDL Astronomy User’s Library 2 . We also fit 
a two-dimensional Gaussian to the stellar PRF (Agol et al. 2010). 
We find that determining the centroid position using GCNTRD re¬ 
sults in a smaller dispersion of the fitted residuals by 10 to 15% 

1 Spitzer Heritage Archive: http://sha.ipac.caltech.edu 
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across our dataset, in agreement with other warm Spitzer analy¬ 
ses (Beerer et al. 2011). We then perform aperture photometry for 
each dataset using a modified version of the APER procedure using 
aperture sizes ranging from 2.2 to 4.4 pixels in 0.2 pixels intervals. 
The background apertures are located 10 to 14 pixels from the pre¬ 
determined centroid position. We also measure the PRF full width 
at half maximum along the X and Y axes. We use a moving av¬ 
erage based on forty consecutive frames to discard datapoints that 
are discrepant by more than 5a in background level, X-Y position 
and FWHM. In average, we find that 0.5% of the datapoints ac¬ 
quired without the PCRS peak-up mode are discarded, while it is 
only 0.06% when the PCRS peak-up mode is enabled. The result¬ 
ing time-series are finally binned per 30 s to speed up the anal¬ 
ysis. The product of the data reduction consists of a photometric 
datafile for each aperture (12 in total per dataset). We then use the 
out-of-eclipse portion of the light curve to measure the RMS and 
level of correlated noise (see Sect. 2.3.2) to determine the aperture 
that minimises both quantities. In the cases of AORs 48070144, 
48073472 and 48073728 we found that two to three aperture sizes 
produced similarly small RMS and beta factors. We selected the 
one with the smallest RMS. For these three AORs, repeating the 
analysis using the other aperture sizes does not change the retrieved 
occultation depths. We show the retained aperture size as well as 
the corresponding RMS and p factor for each dataset in Table 1. 
For all datasets we discard the photometry corresponding to the 
first 30 minutes as it usually exhibits a strong time-dependent ramp 
(Demory et al. 2011). This procedure is common practice for warm 
Spitzer data (e.g. Deming et al. 2011). 

2.3 Photometric Analysis 

We split the analysis in three steps. In a first step we focus on the 
modelling of the intra-pixel sensitivity and the photometric preci¬ 
sion. The second step focuses on an improved determination of the 
system parameters using six transits and eight occultations. The 
third step details the analysis focusing on the variations in transit 
and occultation depths. 

2.3.1 Intra-pixel sensitivity correction 

Warm Spitzer data are affected by a strong intra-pixel sensitiv¬ 
ity that translates into a dependence of the stellar centroid posi¬ 
tion on the pixel with the measured stellar flux at the 2% level 
(Charbonneau et al. 2005). A usual method is to detrend the mea¬ 
sured flux from the X/Y position of the centroid, using quadratic 
order polynomial functions. In the past years, several stud¬ 
ies (Ballard et al. 2010; Stevenson et al. 2012; Lewis et al. 2013; 
Deming et al. 2014) have employed novel techniques that improved 
the mitigation of the intra-pixel sensitivity in Spitzer! IRAC 3.6 and 
4.5/im channels. 

It has been shown that in some cases, the technique used to 
correct the intra-pixel sensitivity could have a significant impact on 
the retrieved system parameters in general and the eclipse depth in 
particular (Deming et al. 2014). In the present work, we use a mod¬ 
ified implementation of the BLISS (BiLinearly-Interpolated Sub¬ 
pixel Sensitivity) method (Stevenson et al. 2012). 

The BLISS algorithm uses a bilinear interpolation of the mea¬ 
sured fluxes to build a pixel-sensitivity map. The data are thus self- 
calibrated. We include this algorithm in the Markov Chain Monte 
Carlo (MCMC) implementation already presented in the litera¬ 
ture (Gillon et al. 2012). The improvement brought by any pixel¬ 
mapping technique such as BLISS requires that the stellar centroid 
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Date [UT] 

Program ID 

AOR# 

AOR duration [h] 

Eclipse type 

Aperture [pix] 

RMS/30s [ppm] 

ft r 

2011-01-06 

60027 

39524608 

4.9 

Transit 

2.6 

359 

1.75 

2011-06-20 

60027 

42000384 

5.9 

Transit 

3.2 

367 

1.45 

2012-01-18 

80231 

43981056 

5.9 

Occultation 

2.8 

364 

1.54 

2012-01-21 

80231 

43981312 

5.9 

Occultation 

3.2 

329 

1.04 

2012-01-23 

80231 

43981568 

5.9 

Occultation 

3.2 

337 

1.00 

2012-01-31 

80231 

43981824 

5.9 

Occultation 

3.2 

362 

1.27 

2013-06-15 

90208 

48070144 

5.9 

Occultation 

2.6 

323 

1.12 

2013-06-21 

90208 

48070656 

8.8 

Transit 

2.8 

365 

1.88 

2013-07-03 

90208 

48072448 

8.8 

Transit 

3.2 

350 

1.16 

2013-07-08 

90208 

48072704 

8.8 

Transit 

2.6 

378 

1.37 

2013-07-11 

90208 

48072960 

8.8 

Transit 

2.6 

400 

1.85 

2013-06-18 

90208 

48073216 

8.8 

Occultation 

3.0 

325 

1.59 

2013-06-29 

90208 

48073472 

8.8 

Occultation 

3.0 

357 

1.33 

2013-07-15 

90208 

48073728 

8.8 

Occultation 

3.4 

381 

1.25 


Table 1. 55 Cnc e Spitzer dataset. Astronomical Observation Request (AOR) properties for 2011-2013 Spitzer/IRAC 4.5-pm data used in the present study. 


remains in a relatively confined area on the detector, which warrants 
an efficient sampling of the X/Y region, thus an accurate pixel map. 
In our implementation of the method, we build a sub-pixel mesh 
made of n 2 grid points, evenly distributed along the x and y axes. 
The BLISS algorithm is applied at each step of the MCMC fit. The 
number of grid points is determined at the beginning of the MCMC 
by ensuring that at least 5 valid photometric measurements are lo¬ 
cated in each mesh box. Similar to Lanotte et al. (2014), we find 
that the PRF’s FWHM along the x and y axes evolve with time and 
allows further improvement on the systematics correction. We thus 
combine the BLISS algorithm to a linear function of the FWHM 
along each axis. We find that including a model of the FWHM de¬ 
creases the Bayesian Information Criterion (BIC, Schwarz 1978) 
for all datasets. Including a linear dependence with time does not 
improve the fit. We find no correlation between the photometric 
time-series and the background level. We finally perform the same 
intra-pixel sensitivity correction using a variable aperture in the re¬ 
duction and including the noise pixel (Lewis et al. 2013) parameter 
as an additional term in the parametric detrending function. As no¬ 
ticed in previous studies (Lewis et al. 2013), we find this approach 
slightly increases the BIC in the 4.5pm channel. We repeat this 
part of the analysis with lightcurves using apertures sizes smaller 
and larger by 0.2 pixels from the one minimising the photometric 
RMS. We subtract the lightcurves using the altered aperture sizes 
to the optimal one and find point-to-point variations embedded in 
the l-cr individual error bars. We further check the incidence of 
binning on the efficiency of the intra-pixel sensitivity correction 
by analysing the 2012-01-18 occultation (obtained with the PCRS 
mode enabled) with binnings ranging from 0.64s (one integrated 
subarray data cube), 6s, 12s and 30s. We find the following depths: 
106±55, 91±57, 111±60 and 87±53 ppm respectively. We thus 
conclude that using binnings up to 30s does not affect the occulta¬ 
tion depth measurement. 

2.3.2 Noise budget and photometric precision 

The data reduction described above computes the theoretical er¬ 
ror for each photometric datapoint. In order to derive accurate un¬ 
certainties on the parameters, we evaluate the level of correlated 
noise in the data following a time-averaging technique detailed in 
Pont et al. (2006). We compute within the MCMC fit a correction 
factor /3 based on the standard deviation of the binned residuals for 
each light curve using different time-bins. We keep the largest cor¬ 
rection factor ft found for each dataset and multiply the theoretical 


errors by this factor. We show on Table 1 the photometric RMS and 
/? for each dataset. We also conduct a residual permutation boot¬ 
strap analysis to obtain an additional estimation of the residual cor¬ 
related noise. We use for this purpose the lightcurve corrected from 
the systematic effects described in Sect. 2.3.1. The resulting pa¬ 
rameters are in agreement with the ones derived from the MCMC 
analysis (Table 2), while their error bars are significantly smaller. 
This result indicates that the error budget is dominated by the base¬ 
line model uncertainties and not by the residual correlated noise. 


2.3.3 Refinement of the system parameters 

In the second step of our analysis, we use as input data to the 
MCMC all six transits and eight occultations. The goal of this part 
of the analysis is to refine the published parameters for the orbital 
and physical parameters of 55 Cnc e. The jump parameters are the 
planetary-to-star radius ratio Rp/R*, the planetary impact param¬ 
eter b, orbital period P, transit duration W and centre To, and the 
occultation depth S occ . In a first MCMC run, we keep the eccen¬ 
tricity fixed to zero as previous studies using radial-velocity data, 
which provides better constraints on yfe sincu, found eccentricity 
values consistent with 0 (Nelson et al. 2014). We also add the limb- 
darkening linear combinations ci = 2«i + U 2 and C 2 = ui — 2 u 2 
as jump parameters, where u\ and U 2 are the quadratic coefficients 
drawn from the theoretical tables of Claret & Bloemen (2011) us¬ 
ing published stellar parameters (von Braun et al. 2011). We im¬ 
pose Gaussian priors on the limb-darkening coefficients and the 
radial-velocity semi amplitude K (Nelson et al. 2014). We exe¬ 
cute two Markov chains of 100,000 steps each and assess their 
efficient mixing and convergence using the Gelman-Rubin statis¬ 
tic (Gelman & Rubin 1992) by ensuring r < 1.01. Results for 
this MCMC fit are shown on Table 2. In a second MCMC run, 
we explore what limits can be placed on the orbital eccentricity 
e and argument of periastron w. We perform a fit with y^ecosw 
and y/esinu; as additional free parameters and find y^ecosu; = 
0.01 !qo 4 and y^esinu; = — 0.05j(o'47, corresponding to a 3- 
sigma upper limit on the eccentricity of e < 0.19. We do not detect 
variations in the transit or occultation timings. 


2.3.4 Transit depth variations 

In this third part of the analysis, we assess whether transit depths 
vary across our 2011-2013 datasets. This MCMC fit includes the 
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Planet/star area ratio R p /R s 
b = a cos i/R* [.R*] 

T 0 - 2,450,000 [BJDtdb] 
Orbital semi-major axis a [AU] 
Orbital inclination z [deg] 

Mean density p p [g cm -3 ] 
Surface gravity log g p [cgs] 
Mass M p [M 0 ] 

Radius R p [R®] 

&OCC (2012) [ppm] 

&OCC (2013) [ppm] 


0 niQ7+0-0007 
U.U18/-0.0007 

0 36” 1 ” 0 ' 07 
u -' 3D -0.09 

5733.008l°;oQ2 

0 01544’*"°' 00009 
u.u±Of*ffc_ 0 .00009 


83 


+2 

-1 

q+0.8 


3^ 

°‘ -0.7 
q qq+0.04 
O.OO-Q 04 

8.08±g;li 

1 09 + O.O 8 
0.08 

47 ±21 
176 ± 28 


Table 2. 55 Cnee system parameters. Results from the MCMC combined fit. Values indicated are the median of the posterior distributions and the 1-cr 
corresponding credible intervals. 


six transit lightcurves only, each of them having an independent 
planet-star radius ratio Rp / R+. To increase the convergence effi¬ 
ciency, we use Gaussian priors on b, P, W and To using the pos¬ 
terior distribution derived from the joint fit in the previous step. 
Results for this MCMC fit are shown on Table 3. The corrected 
lightcurves for each transit are shown on Fig. 1 and 2 while the 
corresponding correlated noise behaviour is shown on Fig. 3 and 4. 

We find that 55 Cnc e’s transit depths at 4.5/rm change by 25% 
in comparison to the mean value (365 ppm) between 2011 and 
2013. Fitting all transits with a straight line yields a reduced chi- 
squared of Xr = 1.2, corresponding to an overall variability sig¬ 
nificance of la, which we consider being undetected. While most 
transit depth values are consistent at the 1 -a level, only the first 
2011 transit appears as an outlier. Converting our transit depths 
to planetary sizes, we find a minimum planet radius of 1.75±0.13 
Earth radii, a maximum radius of 2.25 ±0.17 Earth radii and an 
average value of 1.92±0.08 Earth radii. This value is significantly 
smaller than the published Spzfzer+MOST combined estimate of 
2.17 ±0.10 Earth radii (Gillonetal. 2012). Alternatively, these 
marginal transit depth variations could be due to steep variations 
of the stellar flux, which is however unlikely for this star (see 
Sect. 2.5). 


2.3.5 Occultation depth variations 

In the following we assess the extent of variability in the planetary 
thermal emission. We first perform an MCMC fit that includes all 
eight occultations observed in 2012 and 2013, each of them having 
an individual depth S occ . All other parameters derived in Sect. 2.3.3 
are included as Gaussian priors in this MCMC fit. Results for this 
MCMC fit are indicated in Table 4. The corrected lightcurves for 
each occultation are shown on Fig. 5 and 6 while the corresponding 
correlated noise behaviour is shown on Fig. 7 and 8. Fitting all eight 
occultations with a straight line yields a mean occultation depth of 
83± 14 ppm and a corresponding Xr = 2.6. 

In a second MCMC fit, we independently fit the occultations 
by season to investigate a difference in depth between 2012 and 
2013. We find that fitting both the 2012 and 2013 occultation depth 
values with a straight line yields Xr = 13.6, meaning that constant 
occultation depths fits the data poorly at the 3.7 -a level. We find 
that the occultation depths as measured with Spitzer at 4.5 pm vary 
by a factor of 3.7j[^'g between 2012 and 2013, from 47±21 to 
176±28 ppm. The resulting lightcurves for each season are shown 
on Fig. 9 and their posterior distributions on Fig. 10. To interpret 
the variability in the occultation depth, we use an observed infrared 
spectrum of 55 Cnc (Crossfield 2012) to estimate the correspond¬ 


ing change in the planet’s thermal emission. Assuming a mean 
radius of 1.92 Earth radii, we find that between 2012 and 2013, 
the planet dayside brightness temperature at 4.5/im changes from 

1365125® K t0 2528^229 K - 

The time coverage of our combined 2011-2013 dataset 
(Fig. 11 ) is too scarce to put firm constraints on the variability 
period of the emission and apparent size of 55 Cnc e, motivating 
further monitoring of this planet. 


2.4 Discrepancy with the Demory et al. 2012 analysis 

We investigate the discrepancy between our updated 2012 oc¬ 
cultation depth, 47±21 ppm, and our previous one published in 
Demory et al. (2012), 131±28 ppm, assuming that the different 
models that we used in the two analyses are at the origin of the 
discrepancy. Indeed, in our 2012 work, we had not yet noticed the 
strong dependency of the measured fluxes on the measured widths 
of the PRF in the x and/or y directions. The improvement brought 
by the addition of the PRF FWHM in the baseline model is de¬ 
scribed in detail in Lanotte et al. (2014) for the case of GJ 436. In 
2012, our modelling of the Spitzer systematics was based uniquely 
on polynomials of the measured x and y position of the PRF centre. 
Furthermore, to reach a decent level of correlated noise in the resid¬ 
uals of our fits, we had to add to our models quadratic functions of 
time and, for three occultations out of four, quadratic functions of 
the logarithms of time (Demory et al. 2012, Table 1). The resulting 
complex time-dependent models could have easily “twisted” the 
light curves to improve the fit, affecting the measured occultation 
depth by a systematic error. To test this hypothesis, we perform 
several MCMC analyses of the 2012 Spitzer light curves (as ob¬ 
tained with a photometric aperture of 3 pixels). For all these anal¬ 
yses, the transit parameters (depth, impact parameter, ephemeris) 
are kept under the control of Gaussian priors based on the result of 
our global analysis presented in Sect. 2.3.3, and assuming a circular 
orbit. 

In our first analysis (Tl), we assume a common depth for all 
four 2012 occultations, and the same baseline models than in the 
present work (Sect. 2.3.1) based on functions of the PRF centre 
position and widths. In our second test analysis (T2), we assume 
again a common depth, and we use as baseline models the complex 
functions used in Demory et al. (2012, Table 1). Our third (T3) and 
fourth (T4) test analyses are similar, respectively, to Tl and T2, 
except that we allow each occultation to have a different depth. 

The results of these tests are summarised in Table 5. By 
comparing Tl and T2, it can be seen that the models adopted in 
Demory et al. (2012) result in a larger and more uncertain occul- 
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Date obs 

Transit depth [ppm] 

2011-01-06 

484±74 

2011-06-20 

2873:52 

2013-06-21 

325±76 

2013-07-03 

365±43 

2013-07-08 

406±78 

2013-07-11 

4333:90 


Table 3. 55 Cnee Spitzer transit depths. IRAC 4.5/jm transit depths derived from the MCMC simulation detailed in Sect. 2.3.4 



Time [BJD-2450000] 


Figure 1. 2011 transit lightcurves. Individual analyses of the Spitzer/IRAC 4.5-//m 201 1 transit data detrended and normalised following the technique 
described in the text. Black filled circles are data binned per 15 minutes. The best-fit model resulting from the MCMC fit is superimposed in red for each 
eclipse. 


tation depth, to significantly larger scatters of the residuals, and to 
much larger values of the BIC. The simpler instrumental models 
used in this present work, based only on measured instrumental 
parameters (PRF position and width) result in a much better fit to 
the data and a reduced and more accurate occultation depth. To 
confirm this conclusion, we use the four 2013 55 Cnc e occultation 
light curves. We discard the real occultations, and inject a fake oc¬ 
cultation signal of 100 ppm, assuming a circular orbit for 55 Cnc e, 
the median of the posteriors of our global analysis for the planetary 
and stellar parameters (Sect. 2.3.3), and a shift of 0.16d for the in¬ 
ferior conjunction timing. We then perform two global analyses of 
the four light curves, keeping again the stellar and planetary param¬ 
eters under the control of Gaussian priors based on the results of our 
global analysis, and taking into account the same 0.16d shift for the 
transit timing. In our first analysis, we use our nominal instrumental 
model based only on low-order (1 or 2) polynomial functions of the 
x- and y-positions of the PRF centre and its x- and y-widths. In our 
second analysis, we omit the PRF widths terms and add quadratic 
functions of time and logarithm of time, as in Demory et al. (2012). 
Table 6 compares the results for both analyses. These results show 
clearly that our updated instrumental model fit the data better and 
result in a more accurate measurement of the occultation depth. 


2.5 Stellar Variability 

55 Cnc benefits from a 11-yr extensive photometric monitoring at 
visible wavelengths obtained with the Automatic Photometric Tele¬ 
scope at Fairborn Observatory. 55 Cnc has been found to be a qui¬ 
escent star showing on rare occasions variability at the 6 milli- 


magnitude level, corresponding to a <1% coverage in star spots 
(Fischer et al. 2008). When detected, this variability has a period¬ 
icity of about 40 days, corresponding to the rotational period of the 
star. The variation in stellar surface brightness along the planet path 
on the stellar disk is directly proportional to the variation in tran¬ 
sit depth (Mandel & Agol 2002; Agol et al. 2010). Conservatively 
assuming an average stellar surface brightness variation of 1% at 
4.5/im cannot be causing the transit and occultation variability of 
20 and ~300% respectively observed for 55 Cnc e. We assume that 
the assessment of variability made in the past on the 11 -yr baseline 
remains valid for our observations. We do not have optical data that 
are contemporaneous with our Spitzer data. The extensive photo¬ 
metric monitoring of the host star, underlining the quiescent nature 
of 55 Cnc, confirms that the source of the observed variability orig¬ 
inates from the planet or its immediate environment. 


3 DISCUSSION 

The observed variability in thermal emission and transit depths of 
the planet seemingly requires an episodic source of 4.5-pm opacity 
in excess of the steady state atmospheric opacity. Thermal emission 
from the dayside atmosphere varies by ~300%, detected at over 4- 
a, leading to temperatures varying between ~ 1300-2800 K assum¬ 
ing a stable photosphere in the 4.5/rm bandpass. Similarly, the am¬ 
plitude of the variation in the planetary radius in our 2013 observa¬ 
tions at 4.5 pm ranges between 0.2-0.8 Earth radii, corresponding 
to a radius enhancement of 11-49 %, albeit detected only at a 1-2 
a level. The amplitudes of these variations require excess absorp¬ 
tion equivalent to ten(s) of scale heights of optically thick material 
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Time [BJD-2450000] 


Figure 2. 2013 transit lightcurves. Individual analyses of the Spitzer/IRAC 4.5-//m 2013 transit data detrended and normalised following the technique 
described in the text. Black filled circles are data binned per 15 minutes. The best-fit model resulting from the MCMC fit is superimposed in red for each 
eclipse. 



Figure 3. 2011 transit - RMS vs. binning. Black filled circles indicate the photometric residual RMS for different time bins. Each panel corresponds to each 
individual transit from 2011. The expected decrease in Poisson noise normalised to an individual bin (30s) precision is shown as a red dotted line. 


absorbing in the IRAC 4.5/.tm band. Local changes in the chem¬ 
ical composition or temperature profiles of the observable atmo¬ 
sphere are unlikely to contribute such large and variable opacity on 
the time-scales of our observations which, nevertheless, need to be 
confirmed by hydrodynamical simulations of the super-Earth atmo¬ 
sphere (e.g. Zalucha et al. 2013; Carone et al. 2014; Kataria et al. 
2015). 

Recent studies motivate the question of whether volcanic out- 
gassing from the planetary surface might contribute to the observed 
variability. Given the high equilibrium temperature of the planet 
(> 2000K) the planetary lithosphere is likely weakened, if not 
molten, for most mineral compositions thereby leading to magma 
oceans and potential volcanic activity, particularly on the irradi¬ 
ated day side (Gelman et al. 2011; Elkins-Tanton 2012). Recent 
inferences of disintegrating exo-mercuries in close-in orbits have 
also used volcanic outgassing as a possible means to inject dust 
into the planetary atmosphere before being carried away by ther¬ 
mal winds (e.g. Rappaport et al. 2012, 2014). However, whereas 
extremely small planets (nearly Mercury-size) subject to intense ir¬ 
radiation can undergo substantial mass loss through thermal winds, 
super-Earths are unlikely to undergo such mass-loss due to their 
significantly deeper potential wells (Perez-Becker & Chiang 2013; 
Rappaport et al. 2014). Thus, ejecta from volcanic eruptions on 
even the most irradiated super-Earths such as 55 Cnc e are un¬ 
likely to escape the planet and would instead display plume be¬ 
haviour characteristic to the solar system (Spencer et al. 1997, 
2007; McEwen et al. 1998; Sanchez-Lavega et al. 2015). The ex¬ 


tent and dynamics of the plumes if large enough can cause tempo¬ 
ral variations in the planetary sizes and brightness temperatures and 
hence in the transit and occultation depths. 

Volcanic plumes may be able to explain the peak-to-trough 
variation in thermal emission that we observe from 55 Cnc e. As¬ 
cending plumes could raise the planetary photosphere in the 4.5 
/rm to lower pressures higher up in the atmosphere where the tem¬ 
peratures are lower, thereby leading to lower thermal emission as 
observed in Fig. 11 . Note that the temperature structure is governed 
largely by the incoming visible radiation and visible opacity in the 
atmosphere, whereas the infrared photosphere is governed primar¬ 
ily by the infrared opacity. Therefore, for the temperature profile 
to remain cooler in the upper-atmosphere the visible opacity in the 
atmosphere and in the rising plume would need to be weaker com¬ 
pared to the infrared opacity. The detailed hydrodynamics of a po¬ 
tential volcanic plume is outside the scope of the current work, but 
we can nevertheless assess its plausibility using simple arguments. 

The observed IRAC 4.5-pm brightness temperature (Tb) of 
the planet varies widely between T m i n = 1273 jjg 4 gK and T max = 
2816j(3ggK. If we assume that the T m i n and T max observed do 
indeed represent the maximum and minimum temperatures in the 
atmosphere, they provide a constraint on the temperature pro¬ 
file of the atmosphere. For extremely irradiated atmospheres the 
temperature profile approaches an isotherm in the two extremes 
of optical depth; i.e. a skin temperature at low optical depths 
and the diffusion approximation at high optical depths (see e.g. 
Madhusudhan & Seager 2009; Spiegel & Burrows 2010; Guillot 
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Figure 4. 2013 transit - RMS vs. binning. Black filled circles indicate the photometric residual RMS for different time bins. Each panel corresponds to each 
individual transit from 2013. The expected decrease in Poisson noise normalised to an individual bin (30s) precision is shown as a red dotted line. 


Date obs 

Socc [ppm] 

T b [K] 

2012-01-18 

87±56 

1767 

2012-01-21 

44±28 

j.oo-L_ 374 

2012-01-23 

39±25 

1273+ 271 

xz, I o_348 

2012-01-31 

82±45 

!720lt?o 

2013-06-15 

212±46 

2816l 3 g 

2013-06-18 

201 ±64 

2729^522 

2013-06-29 

I69±62 

2472 ttll 

2013-07-15 

101 ±52 

1894^07 


Table 4. 55 Cnee Spitzer occultation depths and brightness temperatures. IRAC 4.5/xm occultation depths derived from the MCMC simulation detailed 
in Sect. 2.3.5. Brightness temperatures are computed using a stellar spectrum, as detailed in the text. 


2010; Heng 2012). Therefore, Ta.mm and Ta.max could represent 
the two extreme isotherms. The pressures corresponding to the two 
limits depend on the atmospheric composition, but generally lie be¬ 
tween Pmax of ~ 1-100 bar and P m i n of ~10“ 3 bar in models of 
irradiated atmospheres. Fig. 12 shows an illustrative model of a P- 
T profile between (T mln , P min ) and (T max , P max ). Note that the 
curvature of the T-P profile and hence the vertical location of the 
intermediate points are arbitrary in this toy model. The more im¬ 
portant matter is that the maximum extent of the plumes can be 
encompassed within ~10 scale heights of the atmosphere, imply¬ 
ing ~ 500 km in an H 2 O-rich atmosphere, and 200 or 300 km 
for a CO 2 or CO atmosphere, respectively. Such plume heights are 
comparable to those of the largest plumes in the solar system as 
discussed below (Spencer et al. 1997; McEwen et al. 1998). How¬ 
ever, depending on their aspect ratios and locations relative to the 
sub-stellar point multiple plumes may be required to influence the 
disk integrated thermal emission observed. 

The chemical composition of the outgassed material is uncon¬ 
strained by our current observations and would depend on the in¬ 
terior composition of the planet. However, any candidate plume 
material is required to have significant opacity in the IRAC 4.5 
^tm bandpass. Gases such as CO 2 , CO, and HCN could pro¬ 
vide such opacity (e.g. Madhusudhan 2012). Using internal struc¬ 
ture models of super-Earths (Madhusudhan et al. 2012), we find 
that the mass and bulk radius of the planet, given by the min¬ 
imum radius observed, are consistent with an Earth-like inte¬ 
rior composition of the planet, i.e. composed of an Iron core 
(30% by mass) overlaid by a silicate mantle and crust. Previ¬ 
ous studies which used a larger radius of the planet required 
a thick water envelope (Demory et al. 2011; Gillon et al. 2012; 
Winn et al. 2011) or a carbon-rich interior (Madhusudhan et al. 


2012), neither of which are now required but cannot be ruled 
out either. However, under the assumption of an Earth-like inte¬ 
rior, which is also the case for Io, the outgassed products could 
comprise of gaseous CO 2 and sulfur compounds (Moses et al. 
2002; Schaefer & Fegley 2005) along with particulate matter com¬ 
posed of silicates (Spencer et al. 1997; McEwen et al. 1998). The 
steady state composition of the atmosphere can be equally diverse 
(Miguel et al. 2011; Schaefer & Fegley 2011 ; Hu & Seager 2014), 
although the stability of a thick ambient atmosphere is uncertain 
(Castan & Menou 2011; Heng & Kopparla 2012). While gaseous 
products could contribute significant absorption in the IR. such 
as the strong CO or CO 2 or HCN bands in the 4.5 /.tm Spitzer 
IRAC bandpass in which our observations were made, gaseous 
and particulate Rayleigh scattering and particulate Mie scattering 
can cause substantial opacity in the visible and IR wavelengths. 
New high-precision high-cadence spectroscopic observations over 
a wide spectral range, e.g. with HST and ground-based facilities at 
present and with JWST in the future, will be necessary to charac¬ 
terise the chemical composition of the plume material. 

On the other hand, if volcanic plumes occur along the day- 
night terminator of the planet they could also cause variations in 
the transit depths. We detect only marginal variation (1-cr) in the 
transit depths. Nevertheless, we briefly assess their potential impli¬ 
cations in this context. The fractional change in the transit-depth 
due to plumes in the vicinity of the day-night terminator is given 
by 8 ~ (YliLi e iH?)/Ao, where Ao is the apparent area of the 
planet disk in quiescent phase, N is the number of plumes, Hi is 
the height of plume i and ti is the aspect ratio of plume i. Con¬ 
sidering all the plumes to be of similar extent and an aspect ra¬ 
tio of unity (e.g. Spencer et al. 1997), the average plume height 
required to explain the variability is given by H ~ Aq5/N. 
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Figure 5. 2012 occultation lightcurves. Individual analyses of the Spitzer/IRAC 4.5-//m 2012 occultation data detrended and normalised following the 
technique described in the text. Black filled circles are data binned per 15 minutes. The best-fit model resulting from the MCMC fit is superimposed in red for 
each eclipse. 




Time [BJD-2450000] 


Figure 6. 2013 occultation lightcurves. Individual analyses of the Spitzer/IRAC 4.5-//m 2013 occultation data detrended and normalised following the 
technique described in the text. Black filled circles are data binned per 15 minutes. The best-fit model resulting from the MCMC fit is superimposed in red for 
each eclipse. 


For plumes distributed over the entire terminator, the disk-averaged 
plume height is given by 7T av ~ A R p , where A R p is the ampli¬ 
tude of the apparent change in the planetary radius, which implies 
.ffav ~ 1300 — 5100 km ~ 0.1 — 0.4 R p . On the other extreme, if 
only one plume is involved then H ~ 1 — 2.2 x 10 4 km ~ 1 — 2 
R p . Thus, depending on the number of plumes, the plume height 
required to explain the maximum variability in transit depths can 
vary between H ~ 0.1 - 2 R p . The required plume heights are gen¬ 
erally larger than those found in the solar system, though not im¬ 
plausible. The highest plumes in the solar system are typically seen 
in bodies with low gravities, such as Io and Mars. The highest vol¬ 
canic plumes on Jupiter’s moon Io (radius, Ri a = 1821 km) have 
heights of ~ 300-500 km or 0.16-0.27 Ri 0 (Spencer et al. 1997; 
McEwen et al. 1998). Recently, the highest plume on Mars was re¬ 
ported at the day-night terminator with a height of ~250 km (0.07 
J^Mars) and a width of 500-1000 km, i.e with aspect ratio over 2:1 
(Sanchez-Lavega et al. 2015). The source of the Martian plume is 
currently unknown though apparently non-volcanic. 

The required plume heights are in principle attainable given 
the planet's orbital parameters. The Hill radius of the planet is given 
by .Rmii = a(M p /3M s ) 4 ^ 3 where a is the orbital separation and 
M p and M s are the masses of the planet and star, respectively. The 
required plume heights are well within the Hill radius of 55 Cnc 
e, J?Hiii ~ 3.5 R p , implying that the ejecta is unlikely to escape 
the system and will fall back to the surface. Furthermore, follow¬ 
ing Rappaport et al. (2012), the sound speed (v 3 ) at the surface of 
the planet is ~1 km/s which is significantly lower than the escape 


velocity ( v esc ) of 24 km/s at the planetary surface thereby prevent¬ 
ing a Parker-type wind to carry the ejecta beyond the Hill sphere. 
On the other hand, the velocity of the ejecta required to sustain 
.ffav ~ 1300 — 5100 km is between 8-17 km/s, which is below the 
escape velocity of the planet (24 km/s) and is significantly lower 
than the maximal ejection velocities of ~300 km/s observed on Io. 

Alternately, the data may also be potentially explained by 
the presence of an azimuthally inhomogeneous circumstellar torus 
made of ions and charged dust particles similar to Io (Kruger et al. 
2003) that could contribute a variable gray opacity along the line- 
of-sight through diffusion of stellar light. While CO gas in the torus 
could provide the required opacity in the 4.5 /.im bandpass, the 
grey behaviour of dust could provide the same opacity both in the 
visible and infrared. Io’s cold torus rotates with the orbital period 
of Jupiter (Belcher 1987). Assuming the same scenario holds for 
55 Cnc e would mean that the torus rotates with the same period 
as the host star (~42 days), which would cause a modulation of 
the transit and occultation depths over similar timescales. For such 
a circumstellar torus, the optical depth of the time-varying torus 
segment occulting the stellar disk governs the observed stellar ra¬ 
dius and hence the transit depth; for the same planetary radius, a 
denser segment causes a deeper transit compared to a lighter seg¬ 
ment. Thus, a denser segment gradually coming into view could 
cause the observed increase in transit depth. Similarly for the oc¬ 
cultation depth, increasing the optical depth of the segment could 
result in a decrease in the occultation depth. Volcanism could pro¬ 
vide a source for the torus’ material akin to Io (Kruger et al. 2003). 
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Figure 7. 2012 occultation - RMS vs. binning. Black filled circles indicate the photometric residual RMS for different time bins. Each panel corresponds to 
each individual occultation from 2012. The expected decrease in Poisson noise normalised to an individual bin (30s) precision is shown as a red dotted line. 
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Figure 8. 2013 occultation - RMS vs. binning. Black filled circles indicate the photometric residual RMS for different time bins. Each panel corresponds to 
each individual occultation from 2013. The expected decrease in Poisson noise normalised to an individual bin (30s) precision is shown as a red dotted line. 


Analysis # 

Occ. depth [ppm] 

BIC 

x 2 

RMS/15 min [ppm] 

T1 

47±21 

2938 

2596.1 

64 

T2 

156±37 

3583 

3068.5 

96 

T3 

101±55 

61±37 

94±33 

152±46 

2983.4 

2592.9 

64 

T4 

220±64 

3618.6 

3056.1 

93 


55 ±38 

132±||° 

'■) 1 —|—55 

ZA '-61 





Table 5. 55 Cnc baseline model test. The table indicates the results of the different baseline model tests conducted to investigate the discrepancy with the 
results of Demory et al. 2012 (D12). Analysis T1 shows the fit results from the present analysis, using the same data from D12 and the baseline model presented 
in this study. T2 uses the same photometry but with the baseline models from D12. T3 and T4 are similar to T1 and T2 respectively but allow each occultation 
depth to be fit independently. 


Model BIC RMS/15 min [ppm] Occ. depth [ppm] 

p(x,y,Wcc,Wy) 2905 64 97±22 

p(x,y,t,r) 3585 118 63±34 

Table 6. 55 Cnee occultation injection test. The table indicate how well our fake 100-ppm deep occultation is recovered by the baseline model used in the 
present study p(x,y,w x ,w y ) and by the baseline model used in D12 p(x,y,t,r). x, y are the PSF positions, t is time, r is the logarithmic ramp and w x ,w y are 
the width of the PSF along the x and y axes respectively. 
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Figure 9. 55Cnee Spitzer/IRAC 4.5/zm 2012 vs. 2013 occultations. The top and bottom lightcurves are the combined occultations for 2012 and 2013 
respectively. Each season results from co-adding four occultations. Data are detrended from instrumental systematics and normalised following the technique 
described in the text. Black datapoints are binned per 15 minutes. Best-fit occultation models are superimposed in red. 


Future observational and theoretical studies are needed to explore 
this alternative. Constraining the nature of the material surround¬ 
ing 55 Cnc e will be a direct probe of the planet surface composi¬ 
tion and may also constrain the formation pathway of hundreds of 
similar planets that have been found so far. 
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Figure 10. 55 Cnee average occultation depths 8 OC c in 2012 and 2013. Posterior probability distribution functions for the 2012 (red) and 2013 (blue) mean 
occultation depths. Both histograms indicate the entire credible region constrained by the Spitzer data using an MCMC fit. 
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Library and the IDL-Coyote Graphics Library. 
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Figure 11. 55 Cnee Spitzer/IRAC 4.5//m occultation and transit depths. Top and bottom panels indicate the occultation and transit depths with time. Left 
panels are for 2011, middle panels for 2012 and right panels for 2013. The grey solid horizontal lines indicate the mean depths and the dotted lines the l-cr 
credible intervals. 
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Figure 12. 55 Cnee arbitrary pressure-temperature profile. Schematic representation of the observed brightness temperatures in the IRAC 4.5/im relative 
to an illustrative temperature profile of the dayside atmosphere of 55 Cnc e. The data show the measured brightness temperatures along with their uncertainties. 
The solid gray curve shows a schematic temperature profile as discussed in Section 3. 
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